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Abstract 

The non-equilibrium tunnel transport processes are considered in a 
square lattice of metallic nanogranules embedded into insulating host. 
Based on a simple model with three possible charging states (±, or 0) of a 
granule and three kinetic processes (creation or recombination of a ± pair, 
and charge translation) between neighbor granules, the mean-field kinetic 
theory is developed, which takes account of the interplay between charg- 
ing energy and temperature and between the applied electric field and the 
Coulomb fields by non-compensated charge density. The resulting charge 
and current distributions are found to depend essentially on the particular 
conditions in a granular layer, namely, in a free area (FA) or in contact 
areas (CA) under macroscopic metallic contacts. Thus, a steady state 
dc transport is only compatible with zero charge density and ohmic re- 
sistivity within FA, but charge accumulation and non-ohmic behavior are 
necessary for conduction over CA. The approximate analytic solutions are 
obtained for characteristic regimes (low or high charge density) of such 
conduction. Also non-stationary processes are considered, displaying a 
peculiar combination of two strongly different relaxation times. The com- 
parison is done with the available transport data on similar experimental 
systems. 

PACS: 73.40.Gk, 73.50.-h, 73.61.-r 



1 Introduction 

Granular (including nanogranular) systems are of a considerable interest for 
modern technology due to their new physical properties, Hke giant magnetore- 
sistance Coulomb blocade or high density magnetic memory |^, im- 

possible for classical materials. The systems of our interest here are composed 
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from nanoscopic (spherical) metallic particles embedded into insulating matrix, 
and various existing experimental techniques result in ordered or disordered, 
mono- or polydispersed, single- or multilayered structures ||^. They are of no 
less interest from the point of view of fundamental study. In particular, trans- 
port phenomena in granular systems reveal certain characteristics which can- 
not be obtained neither in classical flow regime (in metallic, electrolyte, or gas 
discharge conduction) nor in hopping regime (in doped semiconductors or in 
common tunnel junctions). Their specifics is mainly determined by the dras- 
tic difference between the characteristic time of an individual tunneling event 
(~ U/e-p ~ lO^^'^ s) and the interval between such events on the same granule 
('^ e/(jd^) ~ 10^'^ s, at typical current density j ~ 10~^ A/cm^ and granule 
diameter d ~ 5 nm) . Other important moments are the sizeable Coulomb charg- 
ing energy Ec ^ /eeffd (typically ^ 10 meV) and the fact that the tunneling 
rates across the layer may be even several orders of magnitude slower than along 
it. The interplay of all these factors leads to unusual macroscopic effects, thus, 
a peculiar slow relaxation of electric charge was recently discovered in experi- 
ments on tunnel conduction through granular layers and granular films |^,|@|. 
For theoretical description of such processes we use the model of a single layer 
of identical spherical particles located in sites of a simple square lattice, with 
three possible charging states (±, or 0) of a granule and three kinetic processes 
(creation or recombination of a ± pair, and charge translation) between neigh- 
bor granules, and even such simple model reveals a variety of possible dynamical 
and thermodynamical regimes, to be presented below. 

The detailed formulation of the model, its basic parameters, and its mean- 
field continuum version are given in Sec. ||. Next in Sec. ^ we calculate the mean 
values of occupation numbers of each charging state for steady state conditions, 
including the simplest equilibrium situation (no applied fields), in function of 
temperature. The analysis of current density and related kinetic equation in 
the out-of-equilibrium case is developed in Sec. |^, where also its simple solution 
is discussed for the FA part of the system. The most non-trivial regimes are 
found for the CA part, as described in Sec. || for steady state conduction with 
charge accumulation and non-ohmic behavior. At least, some non-stationary 
solutions are considered in Sec. |^, either obtained straightforwardly from the 
time-dependent kinetic equation or estimated from a simplified equivalent cir- 
cuit. The general integration scheme for non-linear differential equation, corre- 
sponding to steady state in CA, and the particular approximations leading to 
its analytic solutions are dropped into Appendix. 

2 Charging states and kinetic processes 

We consider a system of identical spherical metalHc nanogranules of diameter 
d, located in sites of simple square lattice of period a within a layer of thickness 
6 ~ a of insulating host with dielectric constant e (Fig. 1). In the charge transfer 
processes, granules can bear different numbers a of electrons in excess (or deficit) 
to the constant number of positive ions, so that the resulting excess charge ae 
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defines a certain Coulomb charging energy Ec. At low enough temperatures, 
which are implied in what follows, the consideration can be limited, besides 
the ground neutral state a = 0, only to single charged states a = ±1. For 
the latter case, Ec was given in the classic paper by Sheng and Abeles ||^ in 
the form Ec = f{s/d)/{ed), where s is the mean spacing between granules 
and the dimensionless function f{z) was estimated for a 3D granular array as 
f{z) = 1/(1 + l/2z). Otherwise, the mutual dielectric response of 3D insulating 
host (e/i — e) with volume fraction / < 1 of metallic particles (e™ oo) can 
be characterized by an effective value Ecff — £/{l~ /). For the planar lattice 
of granules, this effective constant can be estimated, summing the own energy 
e^/ed of a charged granule at site with the energy of its interaction with 
electric dipolar moments ~ {e/£eff){d/2n)^n, induced in each granule at site 
n = a{ni,n2): 

Ec = 



1 



e2 



Eeffd 



(1) 



Here the constant a = (a^/8) X^n^^o ~ 0.92, and the resulting e^ff = e(l + 
a{d/a)'^). However, Eq. |l| may underestimate considerably the most important 
screening from nearest neighbor granules at d ^ a, and in what follows we 
generally characterize the composite of insulating matrix and metaUic granules 
by a certain Seff = e'^/dEc ^ £. 

Following the approach proposed earlier |^, we classify the microscopic states 
of our system, attributing the charging variable dn with values ±1 or to each 
site n and then considering three types of kinetic processes between two neighbor 
granules n and n + A (Fig. 2): 

1) electron hopping from neutral n to neutral n + A, creating a pair of 
oppositely charged granules: {an = 0, dn+A = 0) ^ {an = — l,crn+A = 1), it is 
only this process that was included in Sheng and Abeles' theory; 

2) hopping of an extra electron (hole) from n to neutral n + A, that is charge 
transfer: (cr„ = l,o-n+A = 0) ^ {an = 0, cr„+A = 1), or {an = -l,crn+A = 

0) {an = 0, CTn+A = -1); 

3) an inverse process to 1), that is recombination of a pair: {an = 1, Cn+A = 

-1) {<Jn = 0, tJn+A = 0). 

Note that all the processes 1) to 3) are conserving the total system charge 
Q = X^n"'"' hence the possibility for charge accumulation or relaxation only 
appears due to current leads. A typical configuration for current-in-plane (CIP) 
tunneling conduction includes two macroscopic metallic electrodes on top of the 
granular layer, forming contact areas (C A) where the current is being distributed 
from the electrodes into granules through an insulating spacer of thickness 6', 
and a free area (FA) where the current propagates between the contacts (Fig. 
3). To begin with, let us consider a simpler case of FA while the specific analysis 
for CA with an account for screening effects by metaUic contacts will be given 
later in Sec. |^. 

The respective transition rates qi for the processes 1) to 3) are determined 
either by the instantaneous charging states of two relevant granules and by the 
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local electric field Fn and temperature T, accordingly to the expressions: 



iI'a - (1 - <^l) (1 - '^n+A) ^ (eFn • A + , (2) 

'?£A=4(l-^n+A)'/'(-eFn-A), 
(3) tJngn+A (gn^n+A ^ 1) . IT A J? ^ 

9„A = ^ "/^ (ecTn+AFn • A - £;c) . 

Here the charging energy is taken positive, Ec, for pair creation, zero for trans- 
port, and negative, —E^ for recombination. The function ip (E) = ojN-pE/ {e^^— 
1) expresses the total probability, at given (3 — l/fcgT, for electron transition 
between granules with Fermi density of states and Fermi levels differing by 
E. The hopping frequency lo — ujae~'^^'' is expressed through the "attempt" fre- 
quency LOa ~ ep/S-, the inverse tunneling length x (typically ~ 10 nm~^), and 
the intergranule spacing s = a — d. Local electric field Fn on nth site consists 
of the external field Fext and the Coulomb field F^°"' due to all other charged 
granules (from FA): 



7\Coul 



A suitable approximation for such a system is achieved with passing from 
discrete- valued functions of discrete argument n = a(ni , 712 ) to their continuous- 
valued mean-field (MF) equivalents a{r) — (a-n)r (mean charge density) and 
p(r) — (cr^)^ (mean charge carrier density), obtained by averaging over a wide 
enough (that is, great compared to the lattice period but small compared to the 
size of entire system or its parts) area around r, the latter being any point in 
the plane. This also implies passing to smooth local field: 

F(r) = Fe.t + / (7(r')^-^rfr'. (4) 

Seffa^ J |r' - r|^ 

and to averaged transition rates g(*^(r, A) = (q^^^ andp(*^(r. A) — (^o-nQn^^ ■ 
These rates just define the temporal evolution of mean densities: 

d(r) = ^ (r, A) - q^^^ (r + A, - A) - 

A 

-Hp(2)(r + A,-A)-p(3)(r,A)] , 
p(r) = J2 [1^'^ ^) + 9^'^ + ^' - 1^'^ ^)+ (6) 

A 

-Hg(2)(r-HA,-A)-g(3)(r,A) 

The set of Eqs. |^|| is complete to provide a continuous description of the 
considered system, once a proper averaging procedure is established. 
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3 Mean-field densities in equifibrium 

We perform the above defined averages in the simplest assumption of no cor- 
relations between different sites: (/nffn') — (/n) (ffn'), n' ^ n, and using the 
evident rules:((T^'^+^)^ = cr(r), {(t'^)^ = p{r). The resulting averaged rates are: 

g(i)(r, A) ^ a"(r)a°(r + A)^ [er(r) -A + E,], (7) 
g(2)(r,A) =(7°(r + A){a+(r)^[-eF(r) • A] +(7-(r)^[eF(r) • A]} , 
p(2)(r, A) = (T"(r + A) {cr+{r)ip [-eF(r) • A] - cr- {r)ip [eF(r) • A]} , 
9(3) (r, A) = {cr"(r)f7+(r + A)^[eF(r) ■A~Ec] + 

+a+{r)a-{r + A)ip[~eF{r) ■ A- Er]} , 
p(3)(r,A) = {<7-{r)a+{r + A)ip[eF{r) -A-Ec]- 
-<T+(r)r7-(r + A)^[-eF(r) • A-S,]} , 

where the mean occupation numbers for each charging state o'^(r) — [p(r) ± 
a(r)]/2 and (T°(r) = 1 — p(r)satisfy the normalization condition: J2i '''^(i") = 1- 
In a similar way to Eq. |5[ we express the vector of average current density 
j (n) at nth site: 

A 

-p(2)(n + A,-A) +p(3)(n,A) 

then its MF extension j(r) is obtained by simple replacing n by r in the ar- 
guments of g^*^ and p^^K Expanding these continuous functions in powers of 
|A| = a, we conclude that Eq. H reduces to usual continuity equation: 

. ^ . . ^ 

ar = V2-jr, 9 

e 

with the 2D nabla: V2 = {d/dx,d/dy). 

We begin the analysis of Eqs. from the simplest situation of thermal 
equifibrium in absence of electric field, F(r) = 0, then Eq. || turns into evident 
identity: (j(r) = 0, that means zero charge density, and Eq. || yields in zero 
current density: j(r) = 0, while Eq. || provides a finite and constant value of 
charge carrier density: 

2 

2 + cxp [pEc/2) 

At low temperatures, l3Ec ^ 1, this value is exponentially small: po « 
2exp{—f3Ec/2), and for high temperatures, PEc ^ 1, it tends as po ~ poo — 
PEc/9 to the limit poo — 2/3, corresponding to equipartition between all three 
fractions cr* (Fig. 4). 



5 



In presence of electric fields F(r) ^ 0, the local equilibrium should be per- 
turbed and the system can generate current and possibly accumulate charge. 
Our further analysis is aimed to demonstrate that the latter process is impos- 
sible in FA, but, supposing for a moment a non-zero charge density o-(r), we 
obtain from Eq. |6|a more general local relation between a{r) and p{r): 

^2/ N ^ [p(r) - po] [p{r) +PQ- 2pop(r)] 
1 - Po 

describing the increase of charge density with going away from equlibrium (Fig. 
5). Now we are in position to pass to the out-of-equilibrium situations, beginning 
from a simpler case of dc current flowing through FA. 



4 Steady state conduction in FA 

In presence of (generally non-uniform) flelds F(r) and densities o'(r), p{r), we 
expand Eq. || up to 1st order terms in |A| — a and obtain the local current 
density as a sum of two contributions, fleld-driven and diffusive: 

j(r) -j/.ew(r)+j,./(r) = G[p(r)]F(r)-ei?[p(r)]V2a(r), (12) 

where the effective conductivity G and diffusion coefficient D are functions of 
local charge carrier density p: 

G (p) = l£!li^ |p^' (0) + 2 (1 - p) ^' {E,) +Pv{E,)\. (13) 
a 

^ ^ {l-p)WiO)-2il p)^iE.)]^ 
a(p-cr^) 

Note that p — a"^ in Eq. |lj is nothing but the mean square dispersion of charge 
density, and it is expressed through p by Eq. |l^. In view of this equation, we can 
also consider G and D as even functions of local charge density cr, and just this 
dependence will be mostly used below. Also G and D depend on temperature, 
through the functions ^p' . The system of Eqs. |lT[jl4| , together with Eq. ^, 
is closed and self-consistent, defining the distributions of (j[v) and p(r) at given 

Substituting Eq. |l^ into Eq. |^ and considering the low temperature limit 
when the diffusion coefficient D is practically independent of cr, we transform 
the kinetic equation to the form: 

a - a^bDA2a = Va [G{a)F] ■ (15) 

e 

where the l.h.s. has an aspect of common 2D diffusion equation and the r.h.s. 
plays the role of source (or drain) function. For the FA case, with the local field 
F defined by Eq. and the only relevant coordinate being x, along the current 
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j(r) = const (Fig. 3), Eq. |l^ under steady state condition {& — 0) leads to a 
(non-linear) integral equation for cr(x): 



L/2 



L/2 



dx" 



-]\dx'. (16) 



Numeric analysis of Eq. |l^ shows that there is no other its solution but the 
trivial one: o-(r) = 0. Hence there is no diffusive contribution to the current, and 
the steady state of FA in out-of-equilibrium conditions has an ohmic resistivity: 



G(po) 



With the parameter values a ~ 5 nm, u) ~ 10^^ s 
meV, T ~ 12 K, used in Eqs. |l|, |l|, we have r ~ 50 cm, which is in a 
reasonable agreement with the measured resistance i? ~ 50 100 Mf2 in a 10 
layer granular sample of ^ 1 x 1 cm area pi. 



1 eV" 



(17) 
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5 Steady state conduction in CA 

The kinetics in CA includes, besides the processes 1) to 3) of Sec. 3, also 4 
additional microscopic processes between nth granule and the electrode (Fig. 6) 
which are just responsible for variations of total charge Q by ±1. The respective 
rates g'*-*, i — 4,. . . ,7, are also dependent on the charging state (cr, p) of the 
relevant granule and, using the same techniques that before for FA, we obtain 
their mean values as: 

q^^^ = {p + a)ij{-U-K), = (p - a) ^ ([/-£;:), (18) 

= (1 - p) ^ (C/ + , q^'^ = (1 - p) {-U + E'J . 

Here the function ip [E) formally differs from if (E) only by changing the pref- 
actor: ui ^ u' — u)ae~'^^^ <^ oj, but the arguments of these functions in Eq. |8| 
include other characteristic energies. Thus, U — eb'Fc is generated by the elec- 
tric field Fc between granule and contact surface. This field is always normal to 
the surface (see Fig. 7) and its value is defined by the local charge density a. At 
least, the charging energy E'^ for a granule under the contact can be somewhat 
reduced (e.g., by ~ 1/2) compared to Ec- Then the kinetic equations in CA 
present a generalization of Eqs. ||,^, as follows: 

a(r) = ^ (r, A) - q(^'> (r + A, - A) - p^^^ (r, A)+ (19) 

A 

(r + A, - A) - (r, A) _ q(4) (r) + gO) (j.) + ^(6) (j.) _ ^(7) (j.)! , 
p(r) = ['^''^ ^) + 9*'^ + ^' - 9^'^ ^)+ 

A 
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(r + A, - A) - q^^) (r, A) - g^^) (r) - g^^) (j.) + ^(6) (r) + q(^) (r) 

The additional terms, by the "normal" processes 4 to 7, are responsible for 
appearance of a normal component of current density: 

L (r) - ^ (r) - <z(5) (j.) „ g(6) (j.) ^ ^(7) (j.)! ^ (20) 

besides the planar component, still given by Eq. ^. But even more important 
difference from the FA case consists in the fact that the Coulomb field here 
is formed by a double layer of charges, those of granules themselves and their 
images in the metallic electrode (Fig. 7). Summing the contibutions from all 
the charged granules and their images (except for the image of nth granule 
itself, already included in the energy E'^), we arrive at the expression for the 
above mentioned field Fc at the contact surface as a local function of the charge 
density cr(r): 

F,(r)=F^°"'(r,z = fe') = -^^^, (21) 

(instead of integral relations, Eqs. in FA). Then, the planar component 
of the field by charged granules Fpi(r) = ¥'-^°'^^{r,z = 0) is determined by the 
above defined normal field Fc through the relation Fpi(r) = 6'V2i^c(r)- The 
density of planar current is jp/(r) = GFp/(r) — eDV2cr(r), accordingly to Eq. 
|l2| , that is both field-driven and diffusive contributions into jp/ are present here 
and both they are proportional to the gradient of cr(r). In the low temperature 
limit, this proportionality is given by: 

jpzr«- ^5Kr + V2ar. 22 

Note that presence of a non-linear function: 

9i.<y)^^pl + <y^-pl-o\ (23) 

defines a non-ohmic conduction in CA. In fact, this function should be given 
by Eq. up to maximum possible charge density \<Jmax \ ~ •\/ 1 — Pq, turning 
zero for \a\ > \amax\ (note that the latter restriction just corresponds to our 
initial limitation to the single charged states, see Sec. ||). In the same limit of 
low temperatures, the normal current density is obtained from Eqs. p^ , pO| as 
jz(r) = GzFc{r) where Gz ~ uo' N-pE'^Eef f / '^tt ■ Finally, the kinetic equation in 
this case is obtained, in analogy with Eq. |^, as: 

a(r) = -^V2-jp,(r)-|-^j,(r). (24) 

This equation permits to describe the steady state conduction as well as various 
time dependent processes. The first important conclusion is that steady state 
conduction in CA turns only possible at non-zero charge density gradient, that 
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is, necessarily involving charge accumulation, in contrast to the above considered 
situation in FA. 

Let us begin from the steady state regime which is simpler, in order to use the 
obtained results later for the analysis of a more involved case when an explicit 
temporal dependence of charge density is included in Eq. |2^. 

We choose the CA geometry in the form of a rectangular stripe of planar 
dimensions L x L', along and across the current respectively. In neglect of 
relatively small effects of current non-uniformity along the lateral boundaries, 
the only relevant coordinate for the problem is longitudinal, x (Fig. 8). In the 
steady state regime, the temporal derivative & in Eq. |2j is zero and the total 
current / = const, defined by the action of external source. Then, using Eq. 

we arrive at a non-linear 2nd order equation for charge density: 

^i9['T{x)]+r}^~k'aix) = 0. (25) 

The parameters in Eq. |2^ are: = {u' E'^) / {ahujk-QTi) and t = T/Ti, where 
T is the actual temperature and Ti = Sne^b' /a^bseff. To define completely its 
solution, the following boundary conditions are imposed: 



da (x) 



dx 

da {x) 



dx 



_ k^b'a {x = 0) 

n T 

(27) 



LeubN-pk^Ti g[a{x ^ L)] + t 

The condition corresponds to the fact that the longitudinal current jx at 
the initial point of CA (the leftmost in Fig. 8) is fully supplied by the normal 
current jy entering from the contact to granular sample, and the condition ^ 
corresponds to current continuity at passing from CA (of length L along the 
current) to FA. 

Let us discuss the solution of Eq. |2^ qualitatively. Generally, to fulfill 
the conditions, Eqs. p6|j2^ , one needs a quite subtle balance to be maintained 
between the charge density and its derivatives at both ends of CA. But the 
situation is radically simplified when the length L is much greater than the 
characteristic decay length for charge and current density: /cL >• 1. In this 
case, the relevant coordinate is ^ = L — so that the boundary condition |6| 
corresponds to ^ = L ^ oo, when both its left and right hand side turn zeros: 



f^(Olf^oo = 0' 



di 



= 0. (28) 



The numeric solution shows that, for any initial (with respect to that is related 
to X = L, Eq. H^) value of charge density a — 0) — ao, there is a unique initial 
value of its derivative da (^) /c?^|j=o = D{ao) which just assures the limits 
while for da (^) /d^|^=o > D{aQ) the asymptotic value diverges as cr (f ^ oo) ^ 
oo, and for da /d^\^=o < D{ao) it diverges as a oo) ~oo. Then, 
using the boundary condition |2^ and taking into account the relation V — Voao 
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following from Eq. ^ with Vo = Aireb' /{seffO?), we conclude that the function 
-D(cto) generates the I-V characteristics for CA: 



hb'D{V/Vo)[g{V/Vo) 



(29) 



where h — eujN-pk-QTi. 

A more detailed analysis of Eq. is presented in Appendix. In particular, 
for the weak current regime (regime I) when ctq ^ ci = y%52po(po + T) <C 1, 
so that g{a) « po + along whole CA, Eq. ^ admits an approximate 

analytic solution: 



cr (^) = doe" 



l + 6(^ 



1-e 



(30) 



with the exponential decay index A = k/y^po 
I- V characteristics for regime I: 



V 



Ro (r) 



T. This results in the explicit 



(31) 



presented by the line 1 in Fig. 9. For V < Vi — aiVo, Eq. 
initial ohmic resistence (dependent on temperature t): 



Ro(t) 



Vn 



hkb'y/po (7 



31 describes the 



(32) 



which turns non-ohmic for F ~ Vi . But at so high voltages another conduction 
regime already applies (called regime II), where ai ^ o-q ^ 1 and one has 
g{a) ~ cr (see Eq. Following the same reasoning as for the regime I, we 

obtain a non-linear /- V characteristics for regime II: 



hkb' 



V 



3/2 



(33) 



presented by the line 2 in Fig. 9. This law is weaker temperature dependent 
than Eq. which relates to the fact that the conductance in regime II is 
mainly due to dynamical accumulation of charge and not to thermic excitation of 
charge carriers. However, such dependence can be more pronounced if multiple 
charging states are engaged, as may be the case in real granular layers with a 
certain statistical distribution of granule sizes present. 

At least, for even stronger currents, when already (Tq ~ 1, the solutions of 
Eq. ^ can be obtained numerically, following the above discussed procedure of 
adjustment of the derivative D{uq) to a given ctq. Such solutions define the line 
3 in Fig. 9, which generally agrees with the lines 1 and 2, but its asymptotics 
turns to /- 

A simple and important exact relation for the total accumulated charge Q 
in CA is obtained from the direct integration of Eq. p5[ 



(34) 
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where the parameter t = l/tp {—E'^) has a role of characteristic relaxation time 
in non-stationary processes, as will be seen below. Assuming its value t ~ 1 s 
(comparable with the experimental observations |^), together with the above 
used values of oj and Ti, we conclude that the characteristic length scale 
for solutions of Eq. ^sjcan reach up to 10'*a 0.1 mm, which is a reasonable 
scale for a charge distribution beneath the contacts. 



6 Time-dependent processes 

Now we can extend our analysis also to a non-stationary situation. Actually, let 
us consider the simplest example of this type, that of interrupting the circuit of 
Fig. 3 at the moment to = 0, after an infinitely long steady state with current 
/. Beginning from this moment, the constant value is not the total current, as 
was the case in Sec. ^, but the total charge in the overall system composed 
by the granular sample itself and electrodes (in fact, this total value is zero 
and it corresponds to the uniform distribution cr(r,oo) = in the asymptotic 
Hmit t oo). We use the same ID geometry and, accordingly to the above 
mentioned absence of charges in FA, limit ourselves only to consideration of CA. 
The initial distribution <j{x), corresponding to the steady state, will degrade at 
t > to, at the expense of recombination of charges accumulated in granules 
beneath the contact with the polarization charges on the contact surface (the 
"image" charges). This process is described by Eq. [2^, subject to the initial 
condition a{x,t = to) = and the boundary conditions: 



(T{L,t) g[a{x,t)]+T da{x,t) 
a [L, t) — h - 



a(0,t) 



T&fc^ dx 
cr(0,t) g[a{x,t)]+T da{x,t) 



I-CV 

(35) 



rbk? dx 



(36) 

x=0 



where V and C are respectively the voltage and capacitance between the con- 
tacts. In similarity with the above considered steady state case, the non- 
stationary Eq. |2j with the boundary conditions, Eqs. |3^^, can be solved 
numerically, in fact presenting an exponential decay with characteristic time 
T (the corresponding plot not shown here) . 
But a more transparent description of dynamical processes in the considered 
system is achieved with using the equivalent circuit (Fig. 10), following from 
the previous analisis. Here R represents the resistence of FA, which is ohmic 
(due to the absence of charges) but strongly temperature dependent. The two 
resistances R' associated to each CA are already non-ohmic, as seen from Eqs. 

due to the accumulated charge Q. The proportionality between Q and 
the current Iqp^ through CA, Eq. together with the simplifying assumption 
of uniform charge distribution over the rectangular CA, permits to express this 



non-ohmic behavior through the function g given by Eq. 23 



^'(^Ca)--77^- (37) 
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Here the resistance scale Rq = aL/(AbL'e^Lp' (0)) can be fixed by a direct com- 
parison with the experimentally measured steady state resistance of the whole 
circuit, and the current scale /q — eLL' / {to^) is about 1 /x.A for L ~ L' ~ 1 
cm and a ~ 5 nm. The respective capacities have the orders of magnitude: 
Ci ~ L (capacity between the electrodes through FA) and C2 ~ LL' jh (capac- 
ity between an electrode and adjacent granular layer), hence Ci ^ C2. Then 
the temporal evolution of the system is given by the equations for voltages V\ 
and V2 on them: 

RC \RC' C ) ' 

The numeric solution to this equivalent system with initial conditions Vi (0) = 
V , V2(0) = {V — /i?i)/2, adequate to the considered discharge regime, is pre- 
sented in Fig. 11. It can be noted that the character of decay is very similar 
to that obtained with the exact kinetic equation, Eq. |2j. Another intersting 
moment is the presence of two "relaxation times" seen in a log plot (lower panel 
of Fig. 11). The initial discharge phase, with the strongest current through CA 
and the lowest resistance R' (by Eq. ^) , has a very fast decay whereas the final 
phase has a much slower relaxation rate, due to the growth of R' . This behav- 
ior also agrees with the experimental observations at switching off the current 
through a granular layer [H, reveahng two strongly different relaxation times, 
the faster one being temperature dependent. 

The relatively simple non-Hnear system, Eqs. |8[ can be also used for treat- 
ing more complicate transient processes. 



7 Appendix 

Let us consider the equation: 

|[5('T)+r]^-fcV = (39) 

with certain boundary conditions o-(O) = do, (^'(O) = ctq, resulting from Eqs. 
For a rather general function g (cr) we can define the function 

/(a)= r gia')da\ (40) 
Jo 

then Eq. ^ presents itself as: 

^-^ = k'aiO, (41) 

where F (^) = /[(t(^)] -|- Tcr(^). Considered irrespectively of ^: 

f{a)+Ta = F, (42) 
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this equation also defines cr as a certain function of F: a — <j{F). Hence it is 
possible to construct the following function: 

Lp{F) = 2 f a{F')dF'. (43) 

"'0 

Now, multiplying Eq. Ellby 2dF/d^, we arrive at the equation: 



-k^% (44) 



2 

d(\d^) dC 

with = ip[F{^)]. Integrating Eq. ^ in ^, we obtain a 1st order separable 
equation for F{^): 

dp 

— = ±fcv^. (45) 

We expect the function F to decrease at going from ^ = into depth of CA, 
hence choose the negative sign on r.h.s. of Eq. 45 and obtain its expHcit solution 
as 

dF' ,^ 

TlWf\ = (46) 

with Fq — f{<Jo) + tctq. Finally, the sought solution for cr(^) = cr[F(^)] results 
from substitution of the function F(^), given implicitly by Eq. into (t{F) 
defined by Eq. Consider some particular realizations of the above scheme. 

For the function g{(j) given by Eq. |2^, we have the explicit integral, Eq. 
in the form: 



F{a)= f{a)+Ta= {t + ^^U^ p2 _ _ j ^ ^ ^2 / VJ^ 

^ (47) 
In the case cr <C po ^ 1 (regime I), Eq. ^ is approximated as: 

3 

F « (po + r) a + (48) 

hence cr(F) corresponds to a real root of the cubic equation, Eq. and in the 
same approximation of regime I it is given by: 

F ( 8F^\ 
a(F)«^— 1-— , (49) 

PQ + T \ erf / 



with (Ti = 4^ Po{po + t)-^- Using this form in Eq. we obtain: 

F^ / _ 4f 
Po + r \ a 



i?2 / - 
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and then substituting into Eq. ES: 



In 



Fn 



l + Jl-{2Fo/a,Y 



F 



Inverting this relation, we define an explicit solution for F{^): 



-AC 



1-e 



(51) 



(52) 



Finally, substituting Eq. |52| into Eq. 
corresponding to Fig. 12. 

For the regime II we have in a similar way 



we arrive at the result of Eq. M 



Fia)^a{T + a/2), 



a (F) w V2F + t2 - r, 



(53) 



^iF)^l[{2F + r^f'~r{3F + r^) 



F^^'-X,^ + 



3t 



with Ai — 2^/''fc/(4\/3), obtaining the charge density distribution (Fig. 13e): 



(54) 



This function seems to turn zero already at ^ = ( y/^o + t — y^)/Ai, but in fact 
the fast parabolic decay by Eq. |5j only extends to ^ ~ such that ~ po, 

and for £,>£,* the decay turns exponential, like Eq. |3^. The /- V characteristics, 
Eq. pa, follows directly from Eq. |j. 
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Fig.l: Square lattice of granules in insulating matrix. 




Fig. 2: Kinetic processes in a granular layer. 
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Fig. 3: CIP conduction geometry. 
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Fig. 4: Equilibrium density po of charge carriers in function of tem- 
perature (solid line). The curve 1 (dashed line) corresponds to low 
temperature assymptotics po ~ 2exTp{—Ec/2kBT), and the curve 
2 (dash-dotted-line) to high temperature assymptotics po ~ Poo — 
Ec/QkBT, converging to the limit poo = |- 




Fig. 5: Relation between local charge density a and charge carrier 
density p (both normalized to the thermal equilibrium value po). 
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Fig. 7: Formation of local electric field by charged granules and their 
images at the surface of metallic electrode (point a) and between the 
granules (point b). 
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Fig. 8: Relations between nth longitudinal (jx) and normal (j^) 
rents in CA. 
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Fig. 9: I-V characteristics for CA, line 1 (short dashes) for regime 
1, Une 2 (dash-and-dot) for regime 2, and line 3 (solid) for numeric 
solution of Eq. 24. 
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Fig. 10: Equivalent circuit for transient regimes of CIP conductance. 
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Fig. 11: Charge relaxation at disconnecting the steady state circuit. 
Log plot (lower panel) reveals two strongly different relaxation times. 
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Fig. 12 : Charge density and current distribution in CA(regime I). 



1 




Fig. 13: Charge density in regime II. A fast decay by Eq.53 is 
changed to a slower exponential law, Eq.29, after density dropping 
below the characteristic value pQ. 
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